# -*- coding: utf-8 -*-

"""
@Datetime: 2019/3/31
@Author: Zhang Yafei
"""
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.datasets.samples_generator import make_blobs


def init_sample():
    """
    第一步：生成测试数据
        1.生成实际中心为centers的测试样本300个，
        2.Xn是包含150个(x,y)点的二维数组
        3.labels_true为其对应的真是类别标签
    """
    # 生成的测试数据的中心点
    centers = [[1, 1], [-1, -1], [1, -1]]
    # 生成数据
    X, label_true = make_blobs(n_samples=150, centers=centers, cluster_std=0.5, random_state=0)
    return X, label_true


class AP(object):
    """ AP聚类 """

    def __init__(self):
        self.Xn = None
        self.Xn_len = None
        self.R = None
        self.A = None
        self.simi_matrix = None
        self.class_cen = None

    def fit(self, data):
        self.Xn = data
        self.Xn_len = len(data)
        # 初始化R、A矩阵
        self.R = np.zeros((self.Xn_len, self.Xn_len))
        self.A = np.zeros((self.Xn_len, self.Xn_len))
        # 计算相似度
        self.cal_simi()
        # 输出聚类中心
        self.class_cen = self.cal_cls_center()

    def cal_simi(self):
        """
        计算相似度矩阵
            这个数据集的相似度矩阵，最终是二维数组
            每个数字与所有数字的相似度列表，即矩阵中的一行
            采用负的欧式距离计算相似度
        :return:
        """
        simi = [[-np.sqrt((m[0] - n[0]) ** 2 + (m[1] - n[1]) ** 2) for n in self.Xn] for m in self.Xn]

        # 设置参考度，即对角线的值，一般为最小值或者中值
        # p = np.min(simi)   ##11个中心
        # p = np.max(simi)  ##14个中心

        p = np.median(simi)  ##5个中心
        for i in range(self.Xn_len):
            simi[i][i] = p

        self.simi_matrix = simi

    def iter_update_R(self, old_r=0, lam=0.5):
        """
        计算吸引度矩阵，即R
               公式1：r(n+1) =s(n)-(s(n)+a(n))-->简化写法，具体参见上图公式
               公式2：r(n+1)=(1-λ)*r(n+1)+λ*r(n)
        迭代更新R矩阵
        :param old_r: 更新前的某个r值
        :param lam: 阻尼系数，用于算法收敛
        :return:
        """
        # 此循环更新R矩阵
        for i in range(self.Xn_len):
            for k in range(self.Xn_len):
                old_r = self.R[i][k]
                if i != k:
                    max1 = self.A[i][0] + self.R[i][0]  ##注意初始值的设置
                    for j in range(self.Xn_len):
                        if j != k:
                            if self.A[i][j] + self.R[i][j] > max1:
                                max1 = self.A[i][j] + self.R[i][j]
                    ##更新后的R[i][k]值
                    self.R[i][k] = self.simi_matrix[i][k] - max1
                    ##带入阻尼系数重新更新
                    self.R[i][k] = (1 - lam) * self.R[i][k] + lam * old_r
                else:
                    max2 = self.simi_matrix[i][0]  ##注意初始值的设置
                    for j in range(self.Xn_len):
                        if j != k:
                            if self.simi_matrix[i][j] > max2:
                                max2 = self.simi_matrix[i][j]
                    ##更新后的R[i][k]值
                    self.R[i][k] = self.simi_matrix[i][k] - max2
                    ##带入阻尼系数重新更新
                    self.R[i][k] = (1 - lam) * self.R[i][k] + lam * old_r
        print("max_r:" + str(np.max(self.R)))

    def iter_update_A(self, old_a=0, lam=0.5):
        """
        迭代更新A矩阵
        :param old_r: 更新前的某个r值
        :param lam: 阻尼系数，用于算法收敛
        :return:
        """
        old_a = 0  ##更新前的某个a值
        lam = 0.5  ##阻尼系数,用于算法收敛
        ##此循环更新A矩阵
        for i in range(self.Xn_len):
            for k in range(self.Xn_len):
                old_a = self.A[i][k]
                if i == k:
                    max3 = self.R[0][k]  ##注意初始值的设置
                    for j in range(self.Xn_len):
                        if j != k:
                            if self.R[j][k] > 0:
                                max3 += self.R[j][k]
                            else:
                                max3 += 0
                    self.A[i][k] = max3
                    # 带入阻尼系数更新A值
                    self.A[i][k] = (1 - lam) * self.A[i][k] + lam * old_a
                else:
                    max4 = self.R[0][k]  # 注意初始值的设置
                    for j in range(self.Xn_len):
                        # 上图公式中的i!=k 的求和部分
                        if j != k and j != i:
                            if self.R[j][k] > 0:
                                max4 += self.R[j][k]
                            else:
                                max4 += 0

                    # 上图公式中的min部分
                    if self.R[k][k] + max4 > 0:
                        self.A[i][k] = 0
                    else:
                        self.A[i][k] = self.R[k][k] + max4

                    # 带入阻尼系数更新A值
                    self.A[i][k] = (1 - lam) * self.A[i][k] + lam * old_a
        print("max_a:" + str(np.max(self.A)))

    def cal_cls_center(self, max_iter=100, curr_iter=0, max_comp=30, curr_comp=0):
        """
        计算聚类中心
            进行聚类，不断迭代直到预设的迭代次数或者判断comp_cnt次后聚类中心不再变化
        :param max_iter: 最大迭代次数
        :param curr_iter: 当前迭代次数
        :param max_comp: 最大比较次数
        :param curr_comp: 当前比较次数
        :return:
        """
        class_cen = []  # 聚类中心列表，存储的是数据点在Xn中的索引
        while True:
            # 计算R矩阵
            self.iter_update_R()
            # 计算A矩阵
            self.iter_update_A()
            # 开始计算聚类中心
            for k in range(self.Xn_len):
                if self.R[k][k] + self.A[k][k] > 0:
                    if k not in class_cen:
                        class_cen.append(k)
                    else:
                        curr_comp += 1
            curr_iter += 1
            print('iteration the {}'.format(curr_iter))
            if curr_iter >= max_iter or curr_comp > max_comp:
                break
        return class_cen

    def c_list(self):
        # 根据聚类中心划分数据
        c_list = []
        for m in self.Xn:
            temp = []
            for j in self.class_cen:
                n = Xn[j]
                d = -np.sqrt((m[0] - n[0]) ** 2 + (m[1] - n[1]) ** 2)
                temp.append(d)
            # 按照是第几个数字作为聚类中心进行分类标识
            c = class_cen[temp.index(np.max(temp))]
            c_list.append(c)
        print(c_list)
        return c_list


def plot(class_cen, X, c_list):
    # 画图
    colors = ['red', 'blue', 'black', 'green', 'yellow']
    plt.figure(figsize=(8, 6))
    plt.xlim([-3, 3])
    plt.ylim([-3, 3])
    for i in range(len(X)):
        d1 = X[i]
        d2 = X[c_list[i]]
        c = class_cen.index(c_list[i])
        plt.plot([d2[0], d1[0]], [d2[1], d1[1]], color=colors[c], linewidth=1)
        # if i == c_list[i] :
        #    plt.scatter(d1[0],d1[1],color=colors[c],linewidth=3)
        # else :
        #    plt.scatter(d1[0],d1[1],color=colors[c],linewidth=1)
    plt.savefig('AP 聚类.png')
    plt.show()


if __name__ == '__main__':
    # 初始化数据
    Xn, labels_true = init_sample()

    ap = AP()
    ap.fit(data=Xn)
    class_cen = ap.class_cen
    print(class_cen)
    # for i in class_cen:
    #    print(str(i)+":"+str(Xn[i]))
    # c_list = ap.c_list()
    # plot(class_cen=class_cen, X=Xn, c_list=c_list)
